Diffusion and conduction in a salt-free colloidal suspension via molecular dynamics 

simulations 



Sorin Baste4!l 

Lawrence Livermore National Laboratory, 7000 East Ave., Livermore, CA 94550 

Molecular dynamics (MD) simulations are used to determine the diffusion coefficients, elec- 
trophoretic mobilities and electrical conductivity of a charged colloidal suspension in the salt-free 
regime as a function of the colloid charge. The behavior of the colloidal particles' diffusion constant 
can be well understood in terms of two coupled effects: counterion 'condensation' and slowdown due 
to the rela:xation effect. We find that the conductivity exhibits a maximum which approximately 
separates the regimes of counterion-dominated and colloid-dominated conduction. We analyze the 
' electrophoretic mobilities and the conductivity in terms of commonly employed assumptions about 

, the role of "free" and "condensed" counterions, and discuss different interpretations of this approach. 
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' Charged colloidal suspensions exhibit a wide range of interesting equilibrium as well as electrokinetic behaviors. On 
: the equilibrium side these include charge renormalization due to strongcounterion screening [H-Qj large Coulombic 
■r^lj- ' effects in sedimentation profiles , highly tunable phase transitions 0, Q , solvation charge asymmetry effects [t!] , et 
caetera. The dynamic behavior may be even richer due to the coupling of hydrodynamic and electrostatic interactions 
[§-11 1 . The motion of charged, mesoscopically-sized particles such as colloids or polyelectrolytes in applied electric 
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fields, generally referred to as electrophoresis, is relevant to numerous applications, from the detection of elementary 



Q charges [12[ to molecular biology )13|1 and nanofluidics [l^]- As a result, the interest in elucidating fundamental aspects 
c/3 of this phenomenon remains high [l^ [l^ , with computer simulations playing a significant role 17 - 1^ . Although the 
focus of such studies is often only the electrophoretic mobility 15, 18|. the behavior of the self-diffusion coefficient and 



electrical conductivity are also important features of these systems [2l| . Electrophoretic mobility and electrical 
^ ' conductivity measurements are in fact routine characterization tools for the charge and number density of colloidal 
particles 0, [2l| , while the determination of the diffusion constant is commonly used to obtain information on their 
^ size [2^ . Since these are the major parameters determining the stability, structure and phase behavior of colloidal 
O . suspensions [l^l, unraveling the mechanisms controlling charge and diffusive transport and their exact relation 
to colloidal charge, size and density is crucial for understanding the properties of charged suspensions and how to 
potentially control them. Numerical simulations are an important avenue for elucidating such connections and new 
methods based on coupling the equations of motion of the colloidal particles with hydrodynamics equations for the 
solvent have recently proved valuable for studying the electrophoretic mobility [ isl , \WL [24| , while Brownian dynamics 
simulations have yielded all transport coefficients ^25]]. Molecular dynamics (MD) simulations are a complementary 
tool for such studies and can directly capture the complex interplay of short and long range interactions as well as 
CO ' solvent structure and hydrodynamics which are relevant for transport processes in colloids, but have been previously 
. limited to single colloidal particle systems |17i, (2^, ^27,]. In this letter we present MD simulation results for the 
electrical conductivity, self-diffusion coefficients and electrical mobilities of a charged colloidal suspension under no- 
, salt conditions. Such suspensions can provide substantial insight for example into surface charge regulation effects 
28| . and they have recently began to be studied more systematically, with results now available for both structural 
231 and electrokinetic properties [2^. 
The system considered here contains solvent particles, colloidal particles with charge —Ze and counterions of 
^ j , opposite unit charge, i.e. it is a salt-less suspension. The inter-particle potentials consist of short range and Coulomb 
' contributions. The short range interactions are based on the inverse-12, 'soft-sphere' potential, 
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u{r) ^e^^j (1) 

, which we truncate and shift at r/do = 2. (We also define u{r) — 00 for r < 0.) They are: 

ucc{r) = u{r ~ 2Rc) (2a) 
ucc{r) = ucs{r) = u{r - Rc) (2b) 
Ucc{r) = Ussir) = Ucsir) = u{r) (2c) 

, where C, c and s denote the colloidal particles, counterions and solvent particles, respectively; Rc is an impenetrable 
colloidal particle core radius. Such potentials have been employed before to model neutral suspensions (sol . [3l| . 
For temperatures ksT ~ e the effective diameters corresponding to these interactions are well approximated by 
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= dcs = do, dc = "^Rc + do, and dc-y = Rc + d^, and satisfy additivity, dc-^ = {dc + d~f)l2; 7 = c, s. The Coulomb 
interactions are given by 

Vafs{r) = (3) 

47reoer r 

where a, P — c, C . Thus, while the solvent size granularity is explicitly accounted for at the microscopic level by 
the short range interactions, the (relative) dielectric constant e^, as is usually the case 2^ 3^, is not. It would be 



difficult to perform large scale simulations with solvent particles carrying explicit dipoles, and we do not expect that 
they would change the results reported here. 

We focused on suspensions with colloid 'volume fraction' 0c — irncdc^ /Q = 0.1, solvent plus counterions 'volume 
fraction' 00 = 7r(?T.c + ns)dQ/6 = 0.35, and with a colloid- solvent 'diameter' ratio dc/do = 10; nsiS = C,c,s) are the 
number densities and satisfy ric = Znc due to charge conservation. We performed MD simulations of this system 
(which may be viewed as a model nanocoUoidal suspension) in the microcanonical (NVE) ensemble for 8 different values 
of the colloidal charge Z : 0, 10, 20, 30, 40, 50, 70, 100 at an average temperature fc^T = e. If we introduce the Bjerrum 
length Xb = e^/47reoerfcsT and Debye screening distance associated with counterions only, Xu — (4TrXBnc)~^ , the 
simulations correspond to Xs/do = 2.32, while Xo/do varied between 4.23 for the Z — 10 simulations and 1.34 for 
the Z = 100 ones. The particle masses were nis = rric = toq and mc/mo — 1000. All simulations were carried out 
with Nc — 50 and Ng + Nc = 175000 in a box with periodic boundaries. The Coulomb interactions were handled 
using the Ewald summation technique with conducting boundary conditions. The time, electrical conductivity and 
electrical mobility units are Iq — doinio/e)^ , ao = 10"'' x e^to/modQ and /io — eto/mo, respectively. 

After equilibration the simulations were run for ~ 5 x 10^ time steps (~ 4000to), accumulating structural as 
well as dynamic information necessary for the calculation of self-diffusion, electrical conductivity and electrophoretic 
mobilities. The self-diffusion coefficients Dg, S = C,c,s were determined using the velocity autocorrelation relation: 
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Ds^ lim - / (V5(0)-V5(r))dr (4) 

t-+oo 3 

with no tail corrections. The conductivity calculation was done by integrating the charge current autocorrelation, 
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: lim — — / (j'(0) •P(r))dr (5a) 
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j'(i) = ^g.v,(0 = eliiit) ~ Zfcit)] (5b) 

i 

fsit) = Y.^u{t) (5c) 

i 

The electrophoretic mobilities of the colloidal particles and counterions can also be calculated in the Green-Kubo 
framework by integrating the charge - diffusion currents correlations [23| : 

Mf = Ihn — i— /*(j''(0) ■3l(.r))dr (6) 
t-)-oo ZNskBi Jo 

Thus a can be written in terms of the mobilities in its canonical form: 

a = e{ncfi^ - ncZfi^) (7) 

It is worth noting that fj,Q is the electrical mobility in suspension and not for an isolated particle, as often studied 
theoretically and measured under very dilute conditions. 

The coUoid-counterion pair correlation functions - Fig. 1 - show strong counterion stratification at the surface of 
the colloidal particles and the formation of what is typically designated as the Stern layer [23|, with a thickness of 
roughly one counterion diameter. This surface 'condensation' effect and the associated chemical equilibrium between 
'condensed' and 'free' counterions leads to colloid-colloid interactions corresponding to a renormalized or effective 
colloid charge Z^ff, which saturates at large ^ We test this scenario by defining Z^ff in the natural way, as the 

charge contained on the average in a sphere centered on a colloidal particle and extending up to the first minimum 
of gccif), i-G- inside the outer boundary of the Stern layer. Other, related definitions are also possible (33| , but they 
yield similar outcomes. The result is plotted in Fig. 1 (inset) and shows that Z^ff so defined exhibits the predicted 
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saturation behavior in the range of 'bare' charges Z covered in these simulations. Being defined in terms of gcc{r), 
Zeff includes the effect of the counterions size. We note that since the system is salt- free and all counterions are 
identical subtle phenomena associated with size asymmetry between the anions and cations such as charge reversal 



and charge amplification [3J] are not expected here. It would be very interesting to study the influence of such effects 
on the electrokinetic properties of charged colloids, since size-asymmetry is common in real situations. We should also 
remark that for such systems the definition of ^e// would likely need to be revised, perhaps by assuming that the 
charge layer extends up to the first minimum of the total number density of ions; this issue requires however further 
careful study. 

The self-diffusion coefficients of the colloidal particles and counterions - Fig. 2 (see also top inset of Fig. 3), are 
both decreasing functions of Z. Their behavior can likely be understood by extending to charged colloidal systems 
arguments usually employed for the self-diffusion coefficients of electrolyte solutions. The classic analysis due to 



Onsager [35[ singles out the most important contribution to self-diffusion to be the relaxation effect, i.e. the drag 
exercised on a moving ion by its lagging, distorted charge atmosphere. This direct, phenomenological treatment is 
particularly suitable for the colloidal particles, which are primarily screened by the small counterions. In this case the 
ionic atmosphere subject to relaxation should reasonably be expected to consist only of the counterions beyond the 
tightly "bound" Stern layer, so an analysis in terms of the effective charge Z^ff is likely appropriate; we outline it 
below following the ideas in (35| . The screening atmosphere of a charged colloidal particle moving with velocity v will 
lag behind a distance ( equal with the one traveled by the particle in the time r that the atmosphere (a charge shell 
of typical size Ad) needs to equilibrate through diffusive redistribution of the screening counterions: C = wr, where 
T oc Xjj/Dc and Dc = ksT/^ (Einstein relation for counterions), with friction coefficient ^ oc rjdc (Stokes relation 
with viscosity ij). This lag or distortion will result in a retarding force between the charged particle and its screening 

will be 



atmosphere, Fr oc (Z^^ rc'^ / X'^er- Since the Stokes drag on the colloidal particle is Fg oc rjdcv, the total drag force 



Ftotai =F,+Fr^ r]dcv[l + !\B{dc/dc)Zlff/\D] (8) 
, where / is a numerical factor. This corresponds to an effective friction coefficient 

£,1 (X r,dc[l + f{dJdc)Zlff\B/\D] (9) 
, which yields, according to the Einstein relation, the self-diffusion coefficient of the charged colloidal particles: 

The Onsager relation for diffusion in electrolytes with equally sized anions and cations (whose derivation requires 
additional assumptions and is similar but not identical with Eq. [TUl) also determines the factor /, e.g. / ~ 0.1 for 
the charge symmetric case [s^. For the present highly asymmetric charged colloidal system we test the relaxation 
effect prediction by assuming / to be a free parameter. The comparison with the MD results, shown in Fig. 2 (inset) 
with / ~ 0.022, indicates that the concept of renormalized charge is suitable for describing the self-diffusive motion 
of charged colloids in conjunction with the relaxation effect. This may provide a convenient avenue for estimating 
the self-diffusion coefficients of such colloids over a wide range of charged states. We note that the relative diffusion 
constant Dc/Dq satisfies the scaling ansatz introduced in [26| for colloidal suspensions in the low-salt regime. Its 
behavior however is different from the standard Einstein relation usually assumed for charged colloids [22]. 

Electrical conductivity (a) measurements are an important means for characterizing the properties of charged 
suspensions [1, [2l|. The MD simulations reveal that a exhibits a maximum as a function of the 'bare' charge Z - 
Fig . 3. Such a behavior has been previously observed for the clectrophoretic mobility of short polyelectrolyte chains 
[lO|, as well as in simulations of a single charged colloidal particle at low volume fraction [27|, |36|. In the present 
simulations the behavior of the counterion and colloid clectrophoretic mobilities - Fig. 3 (bottom - inset), indicates that 
the conductivity maximum roughly separates the regimes of counterion-dominated and colloid-dominated conduction. 
Theories for electrical conduction in classical charged systems such as electrolytes have a long history (37l . but their 
generalization to charged colloids is difficult. Common interpretations of electrical conduction rely on the counterion 
condensation effect and the assumption of "free" and "condensed" counterion populations that contribute differently 
to charge transport [40i]. We examine below to what extent this picture can be reconciled with the observed behavior. 

Using charge conservation the conductivity equation Eq. [7] can be written as 

a = encZifif (11) 
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, where for simplicity we denote (— A*^) by a positive quantity. On the other hand the assumption of "free" and 
"condensed" counterions is usuahy used to write the conductivity as 



a = e7ic^://(Mc/+/^c) (12) 
j^j^ i^uuiiijcinjiio cixc assumed to be "free", with mobihty /i^. 



where now only ncZ*ff counterions are assumed to be "free", with mobility ^ff. This equation holds exactly 



if {Z — ^eff) counterions are essentially "attached" to each colloidal particle, and therefore have the exact same 
mobility with it. In reality however the "condensed" and "free" counterions are not separated populations, and the 
physical, commonly accepted picture is one of chemical equilibrium between them The present simulations are 
consistent with this picture, as they show a continuous exchange between the "condensed" and "free" counterions as 
previously defined based on gccir). We estimate for example that the average exchange rate (per particle) is of the 
order 3.0 x IO^^Iq^. Thus, as opposed to fif and /i^, which are given by Eq. [6l the "free" ions mobility is 
not a very well defined quantity. Moreover, even if some reasonable interpretation for fif^ is adopted, e.g. the single 
particle, high dilution mobility, ^e// '^^^ only be understood as being defined by Eqs. [Tl] and [121 Thus, there is 
no a priori reason why Z*j:j: should be exactly identified with Zg/ / , although there is a reasonable expectation that 
they should behave similarly. We also note that simple assumptions about the relation between fj,f and /i^ such as 
/if = {Z*j:j/Z)fifj: are not compatible with Eq. [TT] unless /i^ is zero, or at least fi^ <C fif- 

We now consider some possible interpretations of the simulation results, using the Nernst-Einstein relation [sij . 
which connects the self-diffusion constants to the mobilities and the electrical conductivity. The simplest application 
of this relation suggests that the conductivity, Eq. [Til should be written as 

C7 - e\.cZiDc + ZDc)/kBT (13) 

, with mobilities 

= eZDc/kBT (14a) 

/if = eDjksT (14b) 

Not surprisingly this is not a good approximation as it deviates rapidly at low Z from the simulation results, and 
predicts a monotonously increasing conductivity - see Fig. 3 (bottom). The main reason for the large discrepancy 
is likely the implicit neglect in the above relation of the correlations between particles and counterions, which is not 
consistent with the existence of a tightly "bound" counterions layer. Such deviations from ideal behavior have also 
been noted in Brownian dynamics simulations of charged suspensions [2^ . 
The Nernst-Einstein relation can also be written starting from Eq. [121 as 

a = e^ncZ:f^{D,f + ZljfDc)/kBT (15) 

, with mobilities 

= eZljPc/ksT (16a) 

^xff = eD.f/ksT (16b) 

This however requires some interpretation on the meaning of Dcf, the diffusion constant of "free" counterions. We 
assume here that Dc{Z = 0) is a reasonable measure for -De/, and identify it with the diffusion constant of the solvent, 
Ds, since the counterions and solvent particles are identical in the limit Z — 0; this quantity is independent of Z, 
and in fact Dc{Z — 10) is already essentially identical with Dg. This assumption is in line with usual conductivity 



modeling based on Eq. [121 [2l|, |40|. We then identify Z*^^ with Z^ff, and plot the results in Fig. 3. This relation 
underestimates the conductivity results, and overestimates the mobility of the colloidal particles. Moreover, the 
qualitative behavior of these quantities is different than the simulations, as they both exhibit plateaus at high Z. We 
note that the observed decrease in mobility at high Z is in agreement with previous simulations, for example those 
reported in ^36*1. 

To make further progress we adopt therefore the following ansatz: we assume that the mobilities of the counterions 
and colloids /ic and fic are proportional with their respective diffusion constants Dc and Dc , and that the conductivity 
relation reduces to the original Ncrst-Einstcin form Eq. 1131 when Z^jf — > Z. The simplest such form is 

a = e'ncZljfiD, + ZljPc)/kBT (17) 
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, with mobilities 

t,^ = e{Z:ff/Z){Dc/kBT) (18a) 
^e{Z:f^/Z)iDjkBT) (18b) 

The nf equation is consistent with the idea that only a fraction of the counterions, Z*j:j:/Z, participate on average in 
electrical conduction as "free" charges. The conductivity equation Eq. [IT] predicts a maximum, but its value is smaller 
than that observed in the simulations if -Z^*// is identified with Z^ff- To get better agreement we adopt a simple 
rescaling, -^e// ~ ^-^^eff- Incidentally this scaling factor is similar with the one previously determined by comparing 
effective colloidal charges from structural and conductivity measurements (2ll |. However, while such a factor may 
signal a difference between static and dynamic effective charges (see also below), we do not regard this particular 
value as carrying fundamental significance. We plot the results in Fig. 4. The agreement is good for the counterions 
mobility and at least qualitative for the conductivity. On the other hand the fi^ relation appears to hold well only 
at high Z. Nevertheless, it is worth studying further the validity of Eq. [T8k in the effective charge saturation regime 
(high Z), since there it would imply (n^)^^ oc Z, which may perhaps be employed for a direct determination of the 



charge [12 1 



The above discussion highlights the inherent limitations of the effective charge representation and "free" and 
"condensed" counterions picture in the context of electrical conduction and electrophoretic mobilities. To the extent 
that it is appropriate to regard Eq. 118b as implicitly defining a dynamic effective ch arg e, it also supports the notion 



that static and dynamic definitions of the effective charge lead to different results [17|, |21[ . To gain further insight 
into why more counterions may dynamically behave as "free", as suggested by the above effective charge rescaling, 
we consider the pair correlation function between colloidal particles and counterions + solvent particles - gcc+siT) - 
Fig. 4, which exhibits strong layering at short distances and on intermediate length scales converges to 1/(1 — ^c), 
corresponding to a higher apparent interstitial fluid density. We define the position (rb) right after the first peak 
of gcc+si^) where this density is reached as a boundary fluid layer, and calculate the effective charge contained in 
this shell. We also determine the charge shell boundary ru corresponding to Z*j:j:, and flnd that it agrees well with 
Tft, particularly at the highest Z's - Fig. 4 (inset). Both are approximately half a counterion diameter smaller than 
ro, the flrst minimum of gcc{'r), which deflnes the Stern layer and Z^ff. This suggests that at high Z, Z*^^ may 
correspond to the boundary fluid layer at the surface of the colloidal particles, which is thinner than the Stern layer. 
Of course such a connection is at this point somewhat speculative and would certainly benefit from additional study, 
for example by analyzing at the MD level the velocity correlations between the colloidal particles and the "condensed" 
counterions under an applied electric field [13, . 

In sum, MD simulations of a charged colloidal suspension in the salt-free regime yield the diffusion constants, mo- 
bilities and conductivity, and help clarify the interplay between charge saturation and slowdown due to the relaxation 
effect in the diffusive motion of colloidal particles. The system exhibits an electrical conductivity maximum as a 
function of the colloid charge, which roughly separates the regimes of counterion-dominated and colloid-dominated 
conduction. We analyze the electrophoretic mobilities and conductivity using commonly employed assumptions about 
the role of the effective charge, "free" and "condensed" contcrions, and discuss different interpretations and some of 
the limitations of this approach. We also find, in agreement with previous observations, that in the effective charge 
saturation regime the effective transported charge appears to be larger than the one determined by the Stern layer, 
and speculate that it corresponds to the boundary fiuid layer at the surface of the colloidal particles. Future simula- 
tions such as the one presented here may help further elucidate this issue. They should also be useful for studing the 
dependence of the diffusion constants, mobilities and conductivity on the colloidal volume fraction, which has only 
recently been adressed using many particle systems [isj . Finally, it may be interesting to study electrokinetic effects 
in charged colloidal dispersions with lower solvent dielectric constants, where clustering effects may play an important 
role (41| . 

This work was performed under the auspices of the U. S. Department of Energy by Lawrence Livermore National 
Laboratory under Contract DE-AC52-07NA27344. 
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FIG. 1: CoUoid-counterion pair correlation function for Z — 10 (dashed line), 50 (dotted line), 100 (solid line). Inset: effective 
colloidal particle charge Zeff as a function of the 'bare' charge Z. 




FIG. 2: Self-difFusion coefficients of counterions (circles) and colloidal particles (diamonds) as a function of colloid charge Z\ 
D* — Dada/ Dodo, Do = d^/to, a — c,C. Inset: Colloidal particle self-diffusion coefficient as a function of the effective charge 
Z^ff (symbols) and relaxation effect relation - Eq. 1101 (solid line). 



9 



40 



20 




40 



20 








30 t/tn 60 








50 Z 



100 



FIG. 3: Top: time-dependent electrical conductivity for Z = 50 (solid line) and Z = 100 (dashed line); inset: time-dependent 
diffusion coefficient (dashed line) and mean-squared displacement (solid line) for colloidal particles with Z = 50. Bottom: 
electrical conductivity - simulations (circles) and Eg. 1131 (dashed line), Eq. [15] with Z*^jj = ^c// (dotted line) and Eg 1 171 
with ^*// = lAZeff (solid line). Inset: electrophoretic mobilities (absolute values) of the colloidal particles (triangles) and 
counterions (diamonds - shifted up 0.1 units); solid lines are Egs. 1181 and dotted line is eZ^ffDc/kBT). 
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FIG. 4: Colloid- (solvent+counterions) pair correlation function gcc+s{r) - dashed line; dotted line corresponds to 1/(1 — 4>c)- 
Inset: first minimum of gccir) - ro (squares), fluid layer boundary - ri (see text) (diamonds), and charge shell boundary 
corresponding to lAZeff - To- (circles) . 



